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ABSTRACT 

The linear hydrodynamic stability of proto-neutron stars (PNSs) is 
considered taking into account dissipative processes such as neutrino transport 
and viscosity. We obtain the general instability criteria which differ essentially 
from the well-known Ledoux criterion used in previous studies. We apply the 
criteria to evolutive models of PNSs that, in general, can be subject to the 
various known regimes such as neutron fingers and convective instabilities. Our 
results indicate that the fingers instability arises in a more extended region of 
the stellar volume and lasts a longer time than expected. 

Subject headings: convection - instabilities - stars: neutron 
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1. Introduction 

According to theoretical models, neutron stars are formed in the aftermath of type II 
and lb supernova explosions, associated to the gravitational collapse of massive stars (8 to 
30M Q ) at the end of their lives. A number of authors have recognized that convection in 
the newly born hot neutron star can play an important role in both enhancing neutrino 
luminosities and increasing the energy deposition efficiency, which might lead to the 
explosion of a supernova (Epstein 1979; Livio, Buchler & Colgate 1980; Colgate & Petschek 
1980; Smarr et al, 1981; Lattimer & Mazurek, 1981; Burrows& Fryxell 1993; Janka & 
Muller 1994). Epstein (1979) pointed out that the negative lepton gradient that naturally 
arises in the outer layers of the PNS after the shock breaks out of the neutrinosphere 
can be convectively unstable. Later on, Colgate & Petschek (1980) suggested that this 
instability could lead to a complete overturn of the core and a strong enhancement of the 
neutrino transport needed for a powerful supernova explosion on a dynamical time scale. 
However, convection in PNSs can be driven not only by the lepton gradient but the entropy 
gradient as well. The suggestion that dissipation of the shock must result in a negative 
entropy gradient due to both neutronization and dissociation (Arnett, 1987) has been 
confirmed in a variety of models by Burrows & Lattimer (1988). As a matter of fact, despite 
different equations of state and differences in the consideration of neutrino transport, 
the development of negative entropy and lepton gradients seems to be common in many 
simulations of supernova models (Hillebrandt 1987, Bruenn & Mezzacappa 1994, Bruenn, 
Mezzacappa & Dineva 1995) and evolutionary models of PNSs (Burrows & Lattimer 1986, 
Keil & Janka 1995, Sumiyoshi, Suzuki & Toki 1995, Pons et al. 1999). 

Concerning numerical simulations, the situation is controversial. Recent hydrodynamic 
simulations by Keil, Janka & Muller (1996) demonstrate that convection can arise in a 
rather extended region of PNSs and may generally last a relatively long time, > Is. Bruenn 
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and Mezzacappa (1994), using the mixing length approximation and Mezzacappa et. al 
(1998) in two dimensional hydrodynamic simulations found only mild convective activity 
in the region near the neutrinosphere. Note that the two-dimensional hydro code used by 
Keil, Janka & Muller is coupled to a radial gray equilibrium diffusion code, that suppresses 
the neutrino transport in the angular direction, essentially underestimating the stabilizing 
effect of neutrino transport. On the contrary, Mezzacappa et. al (1998) code works at 
the limit in which the neutrino transport in the angular direction is fast enough to render 
the neutrino distributions spherically symmetric, therefore overestimating the stabilizing 
effect. The answer is obviously at some point in between, that only a future self-consistent 
multidimensional calculation can determine. 

Helpful insights into the nature and growth rates of fluid instabilities have been 
achieved thanks to semi-analytical investigations (Grossman, Narayan & Arnett, 1993; 
Bruenn & Dineva 1996), and this work is in that line. The presence of convection in 
PNSs is usually argued by the fact that the necessary condition of instability (the Ledoux 
criterion) is fulfilled in a some fraction of the stellar volume. (Epstein 1979, Livio; Buchler 
& Colgate 1980; Keil, Janka & Muller, 1996). However, the Ledoux criterion may have 
no bearing at all on the physics of convection. If dissipative effects caused by viscosity or 
energy and lepton transport are taken into account, the Ledoux condition does not apply 
and has to be substituted by the appropriate criteria. In the present paper, we derived the 
criteria of instability in PNSs employing the diffusion approximation for neutrino transport 
(Imshennik & Nadezhin 1972, Pons et al. 1999) and show that dissipative processes can 
substantially change the picture. 



- 5 - 



2. Basic equations 



Consider the condition of instability in a plane-parallel layer between z = and z = d 
with the gravity g directed in the negative z-direction. We neglect a non- uniformity of g as 
well as general relativistic corrections to hydrodynamic equations. Calculations show that 
initially the instability arises in a surface layer of PNS thus the plane-parallel approximation 
seems to be accurate, at least, during the initial evolutionary stage and, should give 
qualitatively correct results during the late stage. We assume that the characteristic cooling 
time scale of a PNS is much longer than the growth time of instability thus it can be treated 
in a quasi-stationary approximation. Since convective velocities are typically smaller than 
the speed of sound, one can describe instability by making use of the standard Boussinesq 
approximation (see, e.g., Landau & Lifshitz 1959). We consider the linear instability 
when the equations governing small perturbations can be obtained by the linearization of 
hydrodynamic equations. In what follows, small perturbations of hydrodynamic quantities 
will be marked by a subscript "1". The linearized momentum and continuity equations read 



where p and p are the pressure and density, respectively, v is the kinematic viscosity. 
We assume that the matter inside a PNS is in chemical equilibrium thus the density is 
generally a function of the pressure p, temperature T and lepton fraction Y = (n e + n u )/n, 
with n e and n u being the net (particles minus antiparticles) number densities of electrons 
and neutrinos, respectively, and n = n p + n n is the number density of baryons. Since in 
the Boussinesq approximation, the perturbations of pressure are negligible because the 
fluid motions are assumed to be slow and the moving fluid elements are nearly in pressure 
equilibrium with surroundings, the perturbations of density (pi) and entropy per baryon 
(si) can be expressed, in terms of the perturbations of temperature (T±) and lepton fraction 



pvx = -Vpi + gpi + pi/Avi, 



(1) 



V-vi = 0, 



(2) 
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Pl~-p(Pj; + 6Y 1 ), (3) 
si ps m B c v — + o-Fi (4) 

where is the mass of the baryon (we neglect the mass difference between protons 
and neutrons), (5 and 5 are the coefficients of thermal and chemical expansion; 
P = —(d\np/d\nT) pY , 5 = — (d\n p/dY) pT , c p = (T '/m B )(ds / dT) pY is the specific heat at 
constant pressure and a = (ds/dY) pT . 

The above equations should be complemented by the equation driving the evolution of 
chemical composition and heat balance. We employ the equilibrium diffusion approximation 
(EDA) which can be reliable during the early stage of evolution of the PNS, when the mean 
free path of neutrino is short compared to the density and temperature length scales. In 
this approximation, the diffusion and thermal balance equations read 

„f = -V.F, (5) 

where \i is the neutrino chemical potential and F and H are the lepton and heat fluxes, 
respectively. Note that in the equation for the heat balance we have neglected viscous heat 
production since this term is second order in velocity. In the EDA, these fluxes are given by 

F = -a T VT - a v Vr], (7) 

H = -b T VT - b v Vr], (8) 

where 77 = n/ksT is the degeneracy parameter of the neutrino gas, k B being the Boltzmann 
constant; ar, a v and br, b v are the coefficients of diffusion and the thermal conductivities, 
respectively. These kinetic coefficients can easily be related to the coefficients D 2 , D 3 and 
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D4 introduced by Pons et al. (1999), a T = TD 3 , a v = TTD 2 , b T = ksTTD^, b v = k B T 2 TD 3 , 
where T = (k B T) 2 k B /Gir 2 c 2 h 3 . 

For hydrodynamic considerations, it is more convenient to express F and H in terms 
of VT and VY~ rather than VT and Vr/. Then, neglecting the contribution of the pressure 
gradient (barodiffusion), we have 



(XrVT + xyVF) 



(9) 



where 



H = -(frVT + £ y VY"), 



f drj\ ( drj 

XT = a T + a \ — \ } Xy = a „ I — 



£ T = b T + b v \ — 



p,Y V / p,T 



Assuming that unperturbed quantities are homogeneous, equations 
linearized and written in the form 

ATi 



Yi + vi-Vy = A 2 



T 



+ Ay AVI 



where 



and 



AVT 



Ti AVT ATi AV 

j7 - Vl ■ — = Kt-^- + KyAYi 



^^(Vs-aVF) = ) Vp-VT 



s,y 



(10) 



(11) 



can be 



(12) 
(13) 

(14) 



is the superadiabatic temperature gradient that gives the difference between the temperature 
gradient of a fluid with constant entropy and composition and the actual temperature 
gradient, and we have introduced the following characteristic conductivities 
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(15) 
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The system formed by equations (0), @, (0) and flT3|), together with the corresponding 
boundary conditions, completely determines the behaviour of small perturbations. For the 
sake of simplicity, we consider the case when perturbations are vanishing at the boundaries 
z = and z = d. Note that other boundary conditions cannot change the main conclusion 
of our paper qualitatively 



3. The dispersion equation 

The dependence of all perturbations on time and horizontal coordinate can be chosen in 
the form exp(7t — ikx), where 7 is the inverse growth (or decay) timescale of perturbations 
and k is the horizontal wavevector. For such perturbations, we have 

A = £-*. (16) 
The dependence on the vertical coordinate should be obtained from equations (H), (@), 



( [L2f) and (|13l), which can be reduced to one equation of a higher order, say for v\ z . The 
coefficients of this equation are constant in our simplified model, therefore the solution 
for the fundamental mode with "zero boundary conditions" can be taken in the form 
Viz — sm(7rz/d). Then, the dispersion equation for the fundamental mode is 

7 3 + a 2 7 2 + ail + «o = 0, (17) 

where 

a 2 = oj u + u T + ujyi 



a\ = lutujy — coty^yt + Wi/(wr + wy) — (^ + 

5 



9 

.2/ 



a-o = Uv{ujtUy - ujty^yt) - ^ g (^Y - ^u Y t) 

p 



.2 



13 



—U L I Ut — -^ujty 
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In these expressions, we introduced the characteristic frequencies 

U) T = K T Q 2 , ijjy = XyQ 2 , 

uj Y t = ^tQ 2 , u ty = k y Q 2 , (19) 
2 2 (5gk 2 (AVT) 2 2 5gk 2 dY 

W » = vQ > ^ = ^T-^^' ^ = -^"^' 

where (AVT) Z is the ^-component of AVT and Q 2 = (vr/<i) 2 + k 2 . The quantities u u , u>t, 
and ujy are the inverse time scales of dissipation of perturbations due to viscosity, thermal 
conductivity and diffusivity, respectively; ujyt characterizes the rate of diffusion caused by 
the temperature inhomogeneity (thermodiffusion), and uo^y describes the influence of the 
chemical inhomogeneity on the rate of heat conduction; uj g is the frequency (or, in the case 
of instability, the inverse growth time) of the buoyant wave; ujl characterizes the dynamical 
time scale of the processes associated with the lepton gradient. 

Equation ( |TTD describes three essentially different modes which generally exist in a 
chemically inhomogeneous fluid. The condition that at least one of the roots has a positive 
real part (unstable) is equivalent to fulfilling one of the following inequalities (see, e.g., 
Aleksandrov, Kolmogorov & Laurentiev 1963) 

a 2 < 0, a < 0, a\a 2 < a . (20) 

Since z/, kt and Ay are positive defined quantities, the first condition a 2 < will never 
apply, and the discussion will be reduced to the other two conditions. 

In the particular case of chemically homogeneous plasma with a "standard" transport 
(ujty = uyt = 0), we have 

(7 + ujy ) [(7 + w„) (7 + wr) - u# = 0. (21) 

The first root, 71 = — u>y, describes a stable diffusive mode. Two other roots correspond to 
the ordinary buoyancy modes one of which can be unstable if uj 2 — uj v ujt > 0. Then, the 
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necessary condition of instability reduces to the Schwarzschild condition (AVT) Z > 0, or in 
terms of the entropy gradient, ds/dz < 0. 

If dY/dz ^ but dissipative effects are negligible (02 = 0, ao — 0), we have 

7 3 + 017 = (22) 

with a\ = —(ujg + uj\). The first root is degenerate in this case and the two other roots are 
given by ±V — °i- The condition of instability is (—0,1) = oj 2 + u;| > 0, or 

that represents the familiar Ledoux criterion. 



4. Instability criteria in the general case. 

Generally, when dY/dz 7^ and dissipative effects are important, the conditions of 
instability are more complex and depend on a horizontal wavevector of perturbations, k. 
The temperature and lepton gradients required for instability can be quite different for 
perturbations with different k and to obtain the criteria we have to find the minimal values 
of these gradients. From the properties of kinetic coefficients we have KyAy — ky^t > 0, 
therefore we can deduce from conditions (^) that the gradients are minimal for 
perturbations with the wavevector k minimizing the quantity Q 6 /k 2 , which is reached at 
k 2 = (n/d) 2 /2. Then, we have 









U 2 , 


mm 





(24) 



The instability criteria can now be obtained from the last two conditions in (|20|), which 
read, respectively, 
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+ -^-(/GrAy - ^yA T ) < , (25) 
^(5(is + \ Y ) + Pk y ) - (—p-) (f3(v + k t ) + 5\ T ) + 

277T 4 

+ -^(k t + \ Y )[(u + k t )(u + X Y ) - \ t k y ] <0. (26) 

In the case of a negligible diffusivity, the first criterion yields the Rayleigh- Taylor 
condition of instability, 5dY/dz < 0. This instability can be responsible, for instance, for 
the salt fingers phenomena in the terrestrial oceans. Therefore, following Bruenn Sz Dineva 
(1996), we can conventionally call the instability associated to condition ( P5|) neutron 
fingers. We refer to a neutron finger unstable region as a region where condition (^) is 
fulfilled but not condition (p6|). In the case of a "standard" transport with small viscosity 
and diffusivity, the second criterion yields the Schwarzschild criterion of convection, 
although convection in PNS is, in general, quite different from the Schwarzschild convection 
and may arise in oscillatory or nonoscillatory regimes. A region is said to be convectively 
unstable when both conditions (25) and (|26|) are satisfied. If condition ( p6|) is satisfied but 



condition (|25[ ) is not, the system is said to be semiconvectively unstable. 



Despite the fact that convective instabilities in PNS were originally argued by applying 
the Ledoux criterion (see, e.g., Epstein 1979, Livio, Buchler & Colgate 1980), it appears 
that the Ledoux criterion is not the valid criterion for instability if we allow for conduction 
of heat, diffusion of particles and/or viscosity. In this case, the true criteria of instability 
are (1251) and (P§), and only in two limiting cases our criteria reduces to the Ledoux form. 



Generally, the region where any of the conditions (|25| , |26|) are fulfilled can be very 
different to that given by the Ledoux criterion. 
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5. Results and discussion. 

To obtain the different thermodynamical derivatives, diffusion coefficients and 
conductivities appearing in the coefficients of the dispersion equation, we used the results 
from numerical simulations of PNS evolution performed by Pons et al. (1999). The results 
discussed in this paper correspond to the model labelled GM3np, with a baryonic mass 
of 1.6 Mq. The simulation was carried out using a Henyey-like spherically symmetric 
evolution code, coupled to a ID neutrino transport scheme in the flux-limited diffusion 
approximation. Details about the code and the calculation of the diffusion coefficients using 
opacities consistent with the underlying EOS (Reddy, Prakash & Lattimer, 1998) can be 
found in Pons et al. (1999). 

Since calculations of viscosity at high densities are unreliable (compare, for instance, 
van den Horn & van Weert 1981, Goodwin & Pethick 1982, Thompson & Duncan 1993), we 
calculate the unstable region for different values of v. In Figure 1 (upper panel) we plot the 
stability regions of a 1.6 M PNS corresponding to case v — 0. In all the simulations we use 
a value for d given by the pressure scale, d —\ dhip/dr Three different situations can 
clearly be distinguished. Initially, there is a small convectively unstable zone (darkest) near 
the surface that lasts for only a few seconds, surrounded by a neutron finger unstable region 
(intermediate gray tone). The inner part of the star is stable (lighter). After a few seconds, 
the neutron finger instability moves inward, occupying a large portion of the PNS, whereas 
the stable region shrinks. Later on, at t — 12 s, the innermost core becomes convectively 
unstable for a few seconds, while the neutron finger unstable region begins to shrink. By 
about t = 30 s, the PNS is mostly stable, except a small region near the center that is still 
subject to the neutron fingers instability. Instability completely disappears after m 40 s. 

To study the effect of viscosity, Figure 1 (lower panel) shows the time dependence 
of the boundaries of unstable regions for v = r] 2 K T (van den Horn & van Weert 1981). 
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Qualitatively, the unstable regions evolve in a similar fashion to the inviscid case. The 
reason why viscosity hardly influences the neutron finger instability region is that the third 
term on the l.h.s. of (p5|) , which describes dissipation of perturbation due to neutrino 
transport and depends on u, is relatively small compared to the first two terms which are 
independent of v. Therefore, the neutron finger unstable region is approximately the same 
as in Fig.l at any evolutionary stage. On the contrary, convective instability is sensitive to 
the value of viscosity through condition (PB"|). The third term on the l.h.s. of ( |2"6D again 
does not yield too much contribution, however, two other terms depend on viscosity. Due 
to this, the region unstable to convection turns out to be more extended at high viscosity. 
Notice that in the case ?] > 1 the viscosity is much greater than the diffusivity and 
thermal diffusivity but, nevertheless, is not sufficiently large to make the term in the r.h.s. 
important, thus recovering Ledoux criterion from (^). Note, however, that for both models 
the convectively unstable region always lies inside the region unstable to neutron fingers. 
In our simulations we have not found any region unstable to semiconvection. Although this 
result might seem different from Bruenn & Dineva (1996) we should keep in mind that 
the thermodynamic conditions used in this work to describe the interior of the PNS are 
quite different from the conditions used by Bruenn & Dineva (1996) to describe the matter 
just below the neutrinosphere of the collapsed core. Bruenn & Dineva (1996) studied the 
instability in a region, typically with a density of 3xl0 13 g/cm 3 and with a high entropy per 
baryon s = 4, since at such early times (~ 50 — 230 ms after bounce) the excess of entropy 
of the shocked mantle has not been yet radiated away. In this work, however, we focus on 
the long term evolution of a PNS, and most of the region of our interest is at supranuclear 
densities (p > 3 x 10 14 g/cm 3 ) and moderate entropies. 

Since the superadiabatic temperature gradient and the gradient of lepton number play 
a fundamental role in the stability criterion we have shown in Figure 2 graphics of these 
quantities for different times in the evolution. We plot thicker lines for the unstable regions 
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according to the model shown in the lower panel of Figure 1. As we can see from the figures, 
at early and middle times in the evolution, instability region corresponds approximately 
to that with a positive value of the superadiabatic temperature gradient. The effect of 
the negative lepton fraction gradient and/or dissipative processes increase only slightly 
the unstable region. At late time however, the unstable region almost disappears while 
the superadiabatic temperature gradient remains positive in a much wider region. This 
is chiefly caused by the change in the sign of the thermodynamic derivative 5, when the 
lepton fraction drops some critical value. 

In the upper panel of Figure 3 we show the maximum of the real parts of the roots 
of equation (|TT| ) at two different stages in the evolution. A positive value of this quantity 
implies instability and, in that case, the instability growth time is calculated by taking its 
inverse (lower panel). As we can contrast in the lower panel of Figure 3 the characteristic 
growth time of instabilities at the beginning of the evolution (t=l s) is of the order of tenths 
of millisecond. Since this time is much shorter than the evolution time scale, convective 
heat transport might be important in this region. At the middle stages of the evolution 
(t=20 s), two different unstable regions can be distinguished, a first one with a growth time 
of the order of 1 ms (convection) and a second region with associated growth times of the 
order of few tens of millisecond (corresponding to the region unstable to neutron fingers). 
The growth time in both cases is, again, much shorter than the time evolution and, in 
principle, instabilities can grow fast enough to influence the heat transport. 

In summary, we have to emphasize that neutrino transport plays a significant role 
in the stability properties of PNSs. Criteria (25) and (26) obtained in the present paper 
correspond properly to the conditions of PNS where heat transport and diffusion can be 
sufficiently fast. These criteria differ essentially from those used by other authors in this 
context. Our analysis shows that the unstable region inside the PNS can be more extended 
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than predicted, for instance, by the Ledoux criterion alone. Generally, neutron fingers 
instability arises in a more extended region of a stellar volume and lasts a longer time than 
convection. Therefore, this instability seems to be as important as convection for the early 
evolution of a newly born neutron star. Our main conclusion is that, surprisingly, diffusive 
effects allow for the existence of new unstable modes and the role of convection in PNS 
should be carefully reexamined. 
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Figure captions 



Fig. 1.- Evolution of the different unstable zones in a M = 1.6M Q PNS. The darker 
zone corresponds to the convectively unstable regions, the lighter zone to convectively 
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stable regions and the intermediate zone presents the secular instability usually denoted by 
neutron fingers. The upper and lower panels displays the results of the inviscid and viscous 
[y = i] 2 kt) models, respectively. 

Fig. 2.- Superadiabatic temperature gradient (left panels) and lepton fraction gradient 
(right panels) for three different times in the evolution. Thicker lines correspond to the 
unstable region according to Figure 1. 

Fig. 3.- Maximum value of the real parts of the roots of equation (17) at two different 
stages in the evolution (upper panel) and growth time for the unstable regions (lower 
panel). 
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Fig. 3.- 



